An efficient multilevel image thresholding method based on improved heap-based optimizer

Image segmentation is the process of separating pixels of an image into multiple classes, enabling the analysis of objects in the image. Multilevel thresholding (MTH) is a method used to perform this task, and the problem is to obtain an optimal threshold that properly segments each image. Methods such as the Kapur entropy or the Otsu method, which can be used as objective functions to determine the optimal threshold, are efficient in determining the best threshold for bi-level thresholding; however, they are not effective for MTH due to their high computational cost. This paper integrates an efficient method for MTH image segmentation called the heap-based optimizer (HBO) with opposition-based learning termed improved heap-based optimizer (IHBO) to solve the problem of high computational cost for MTH and overcome the weaknesses of the original HBO. The IHBO was proposed to improve the convergence rate and local search efficiency of search agents of the basic HBO, the IHBO is applied to solve the problem of MTH using the Otsu and Kapur methods as objective functions. The performance of the IHBO-based method was evaluated on the CEC’2020 test suite and compared against seven well-known metaheuristic algorithms including the basic HBO, salp swarm algorithm, moth flame optimization, gray wolf optimization, sine cosine algorithm, harmony search optimization, and electromagnetism optimization. The experimental results revealed that the proposed IHBO algorithm outperformed the counterparts in terms of the fitness values as well as other performance indicators, such as the structural similarity index (SSIM), feature similarity index (FSIM), peak signal-to-noise ratio. Therefore, the IHBO algorithm was found to be superior to other segmentation methods for MTH image segmentation.

Segmentation has an important role in the field of image processing 1 . Segmentation is the process of separating an image into two or more homogeneous segments based on the characteristics of the pixels in the image. It is utilized in various scopes, such as industry and medicine 2 , agriculture 3 , and surveillance 4 . Thresholding is one of the most common image segmentation approaches. To define the thresholds, most methods use the histogram of the image 5 , which is vital for determining the probability distribution value of pixels in the image 6 . Thresholding obtains the information of the histogram from an image and determines the best threshold ((th)) for categorizing the pixels into various groups. Image thresholding approaches can be categorized into two types: multi-level and bi-level thresholding. Bi-level thresholding techniques use one threshold to separate an image into two groups, whereas multi-level thresholding (MTH) uses two or more thresholds to separate an image into many groups 1 .
To obtain the best threshold values in MTH segmentation, thresholding techniques can be classified into two approaches: non-parametric and parametric. In parametric techniques, each group of grayscale range should be consistent with a Gaussian distribution. Parametric approaches are dependent on the evaluation of the histogram using mathematical operations. The Gaussian mixture is widespread, where used to define the set of operations that convergent the histogram, and the best thresholds are then selected. Non-parametric approaches employ distinct methods to separate the pixels into homogeneous areas; then, the best threshold is defined using statistical information, such as entropy or variance. The Kapur method 7 and Otsu method 8 are used in this study. The Otsu method selects the best thresholds by the maximization of the variance among groups. In the Kapur method, the threshold value is defined by minimizing the cross entropy between a segmented image and the original image. These methods are efficient for one or two th values of thresholds. However, they have several restrictions; for example, they are very costly in computation, mostly when the number of thresholds increases. Non-parametric techniques have several advantages. Specifically, in terms of computation, these methods are computationally faster than parametric methods, especially when used in optimization problems. Metaheuristic www.nature.com/scientificreports/ with many levels of complexity. The segmentation results are estimated using various assessments, such as the structural similarity (SSIM) index 48 , feature similarity (FSIM) index 49 and peak signal-to-noise ratio (PSNR) 50 . Furthermore, IHBO algorithm was evaluated on the CEC'2020 test suite and compared against seven wellknown metaheuristic algorithms including the basic HBO 29 , SSA 51 , MFO 52 , GWO 53 , SCA 54 , HS 55 , and EMO 27 . The evaluations are executed through various non-parametric and statistical techniques to determine whether the optimal solutions provided by the IHBO are superior.
The main contributions of this paper can be summarized as follows: • An efficient HBO based on OBL called IHBO to overcome the weaknesses of the original HBO is presented.
• Evaluating the effectiveness of IHBO on the CEC'2020 test suite.
• IHBO is proposed to solve the problem of high computational cost for MTH .
• Proving the ability of the IHBO to solve the image segmentation problems using the Kapur's entropy and Otsu's method as fitness function. • Verify the image quality using set of metrics FSIM, PSNR and SSIM to obtain optimal solutions. • Evaluating the performance of the provided method based on the various segmentation degrees to estimate stability of the optimizer and evaluate quality of the segmentation.
The remainder of this paper is organized as follows. "Preliminaries" section describes the materials and methods used in this study, while "The proposed IHBO algorithm" section presents the proposed algorithm. "Environmental and experimental requirements" section illustrates the environmental and experimental requirements, while "Experimental results and discussion" section presents the performance evaluation and experimental results. Finally, conclusions and proposals for future work are provided in "Conclusions and future works" section.

Preliminaries
This section introduces the materials required to implement the proposed segmentation method, as well as the approaches implemented based on the above-mentioned approaches.
Objective functions formulation. The entropy criterion of the Kapur 7 approach and between-class variance of the Otsu 8 approach are widely utilized to determine the optimal threshold value th in image segmentation. Both algorithms were developed for bi-level thresholding techniques. An approach can be readily extended for solving MTH problems.
Otsu method for segmentation. The Otsu method is an automatic and non-parametric technique used to determine the optimal thresholds of an image 8 . This method is based on the maximum variance of the various classes as a criterion to segment the image. The intensity levels L are taken from a grayscale image, and the equation below is used to calculate the probability distribution of the intensity value: where i is a specific intensity level in the range 0 ≤ i ≤ L − 1 and n i is the number of gray level i appearing in the image. The number of pixels in the image is nk and Ph i is the probability distribution of the intensity levels. For the simplest segmentation (bi-level), two classes are represented as The probability distribution for C 1 and C 2 are ω 0 (th) and ω 1 (th) , respectively, as illustrated in (3).
Here i represents a class, and ω i is the occurrence probability, and µ j is the mean of a class. For MTH, these values are obtained as and Kapur entropy. Another non-parametric method used to determine the best threshold value of an image was proposed by Kapur in 7 . The approach determines the best (th) implying the overall entropy to be maximized. For a bi-level scenario, the Kapur target capacity can be determined as where the entropies H 1 and H 2 are computed as follows: In (13), Ph i is the probability distribution of the intensity levels, which is computed by (1), and ω 0 (th) and ω 1 (th) are the probability distributions of classes C 1 and C 2 , respectively. ln(.) represents the natural logarithm. Like the Otsu method, the entropy-based method can be modulated for MTH values. In this case, it is necessary to separate an image into n groups using a similar number of thresholds. The equation below can define the new objective function: where TH = [th 1 , th 2 , . . . th n−1 ] is the vector including MTH. Each entropy is computed separately with its respective th values; thus, (14) is expanded for n entropies as follows: Therefore, the values of probability occurrence (ω c 0 , ω 1 , . . . , ω n−1 ) of n classes can be determined using (10) and the probability distribution Ph i in (1).

Heap-based optimizer (HBO).
The HBO mimics the job responsibilities and descriptions of the employees within a company 29 . Although the job title differ from company to another and from business to another, they are organized in a hierarchy and many of titles are given like corporate hierarchy structure, organizational chart tree, or corporate rank hierarchy (CRH), etc. The collection of methods that outlines how particular activities are directed to realize the goals of an organization and also defines how information flows among levels within the company 56 is called an organizational structure. In this section, we explain the mathematical model of the Heap-based optimizer.
Mathematical modeling of the interaction with immediate boss. The upper levels set the rules and laws for employees within the centralized structure and subordinates follow their immediate boss. By the assumption that each immediate boss is a parent node of its children, thus we can model this behaviour by upgrading the location of each search agent x i with reference to its original node B by using the below equation: (10) ω n−1 (th) = L i=th n +1 Ph i (11) µ n−1 = L i=th n +1 iPh i ω 1 (th n ) .
(12) F kapur (th) = h 1 + h 2 , Ph i ω 1 ln Ph i ω 1 . www.nature.com/scientificreports/ where t is the current iteration, and | | calculates the absolute value. k is the k th component of vector , and it is generated random as following where r is a random number in range [0, 1] . In Eq. (16), the designed parameter is γ , this parameter is computed by the following rule: The current iteration is t, T is the maximum iteration's number, and C is a user defined parameter. while executing the iterations, γ decrease linearly from 2 to 0 and when reach to 0, it will increase again to 2 with iterations.
Modeling the interaction between colleagues mathematically. The employees having the same position are considered to be colleagues. Each employee interact with others to achieve the goals of an organization. By assuming that the nodes at the same level in heap are colleagues and others are search agents x i and they update their position based on the position of others selected colleagues S r , the position of a search agent is calculated as follows: where f is the objective function and calculates the fitness of each search agent. Equation (19) enables the search agents to explore the search space S k r if ( S r ) < f ( x i (t)) and allows to explore the search space x k i otherwise.
Self contribution of an employee. This stage explains the concept of employees self contribution. Modeling of this behavior are executed by retaining the prior position of the employee in the next iteration, as illustrated in below equation: In Eq. (20), the search agent x i does not change its rank for it's kth design parameter in the next iteration. We used this behavior to organize the rate of change of each search agent in population.
Putting it all together. This phase explains how to combine the equations of position updating and modelling in previous subsections in one equation. There are three probabilities of selection that are used to determine equation used in updating position of search agents, this probabilities of selection is used to switch between exploration and exploitation phase. These probabilities is divided into three proportions p 1 , p 2 , and p 3 . The search agent updates its location using Eq. (20) according to the proportion p 1 . The below equation computes the outlines of p 1 .
The current iteration t, T is the maximum number of iterations. The search agent updates its location using Eq. (16) according to the selection of proportion p 2 . The below equation compute the outlines of p 2 .
Finally, the search agent updates its location using Eq. (19) according to the selection of p 3 . The below equation computes the outlines of p 3 .
A general position updating mechanism of HBO is computed as follows: where p 1 , p 2 and p 3 are random numbers inside range [0, 1]. This subsection clarifies that the Eq. (20) improves exploration phase, Eq. (16) improves exploitation phase and convergence, and Eq. (19) allows the search agent to move from the exploration phase to exploitation phase. According to this observations, p 1 is higher initially and decreases linearly over iterations, this decreases the exploration phase and improves exploitation phase www.nature.com/scientificreports/ with iterations. After calculating p 1 , the remainder of the span is splitted into two equal portions, which makes attraction towards the colleague and boss equally probable.
Steps of HBO. This section summarizes the HBO steps and clarifies details about their implementation-related calculations.
• Parameters initialization and definition: At first, all the search agents are randomly initialized in a potential solution space. The minimum and maximum boundaries of the search space are defined by variables (L i , U i ) respectively. The number of the population is (N) and maximum number of iteration (T). The specific parameter C can be calculated from C = ⌊T/25⌋. • Population initialization: The random population P is generated from N search agents, each consisting of D dimensions. The population's representation P is shown as follows: • Heap building: We utilize 3 − ary heap to execute CRH. Although heap is a tree shaped data structure, it can be executed using an array. The below operations are d − ary heap based operations required for the HBO execution.
1. parent (i): By the assumption that the heap is performed as an array, this method receives the node's index then retrieves its parent's index. The formulation of parent's index for a node i is calculated by below equation: where ⌊⌋ indicates the floor function, which retrieves the highest integer less than or equal to a given input. 2. child (i; k): The node can own a maximum of 3 childrens in a 3 − ary heap. Therefore we can say, the manager may not manages more than 3 direct persons. The index of the kth child of a node i is returned by this function. The below equation shows mathematical formulation of this function.
For example,index of the 3nd child of 3nd node is calculated as: 3. depth (i): Assuming the last level depth equals to 0, therefore we can calculate the depth of any node i in constant time through below formula: The ceil function is ⌈⌉ , which retrieves the smallest integer greater than or equal to the input. For example, depth of a node indexed 27 in heap is calculated as: depth(27) = log 3 (81 − 27 + 1) − 1 = ⌈2.6476⌉ = 3 4. colleague (i): Assuming that nodes at the same level of node i are the colleagues of this node. The index of any elected colleague of node i is returned by this step and the index can be calculated by generating any random integer in the range dd depth(i)−1) −1

Heapify_Up (i)
: searching upward in the heap then add node i at its correct place to save the heap property.
Algorithm 1 show the pseudo code of this operation.
Finally, the algorithm to build the heap is shown in Algorithm 2. www.nature.com/scientificreports/ 6. Repeated applications of position updating mechanism: search agents' position is repeatedly updated according to previously explained equations trying to converge on the optimum global. The pseudo code of HBO is shown in Algorithm 3.
Opposition-based learning (OBL). The idea of opposition-based learning (OBL) is applicable strategy of search strategy to avoid stagnancy in candidate solutions. OBL is a novel concept inspired from the opposite relationship between entities 57 . The concept of opposition was presented in 2005 as the first time, which has attracted a many of research efforts in the last decennium. Many of Met-heuristic algorithms use the concept of OBL to develop their performance such as harmony search algorithm 58 , grasshopper optimization 59 , ant colony optimization 60 , artificial bee colony 61 and etc. OBL improve the exploitation phase of a search mechanism. Mostly in meta-heuristic algorithms, convergence occurs quickly when the initial solutions are closer to the optimal location; moreover, late convergence is expected. So that, OBL method produce novel solutions by considering opposite search areas which may prove to be nearer to the best solution. OBL is regraded not only the candidate solutions obtained by a stochastic iteration scheme, but also their ' opposite solutions' located in opposite parts of the search space. The OBL method has been hybridized with many bio-inspired optimization gives shorter expected distances to the best solution compared to randomly sampled solution pairs 62 such as cuckoo optimization algorithm 63 , shuffled complex evolution algorithm 64 , particle swarm optimization 65 , harmony search 66 , chaotic differential evolution algorithm 67 , and shuffled frog algorithm 68 . In optimization problems, the strategy of simultaneously examining a candidate and its opposite solution has the purpose of accelerating the convergence rate towards a globally best solution. According to previous related works, in initialization phase utilize OBL only to improve the convergence and prevent stuck in local optima of HBO, then IHBO is utilized to solve problem of multi-thresholding for image segmentation by use two objective functions called Kapur and Otsu.

The proposed IHBO algorithm
In this paper, the HBO algorithm is enhanced based on the OBL as local search strategy called IHBO to evade the drawbacks of the random population and improve the rate of convergence of the algorithm by developing the variety of its solutions. IHBO uses OBL strategy in the initialization phase to improve the search process as following: Initialization phase. In this phase, the algorithm starts by reading the image, converting it to grayscale, computing the histogram of the selected benchmark images, and then computing the probability distribution by (1). The algorithm initializes the IHBO parameters, which are the population size (N), maximum iteration number (T), boundaries of the search space ( L I , U I ), and number of iterations per cycle (t). Thereafter, the OBL strategy is utilized to calculate the Q i vector-maintaining solution by (28).
Updating phase. This phase provides the best threshold values by evaluating the fitness value of X i and Q i populations. To update the search agents' positions (X), we use the fitness value of the optimal threshold of the Otsu F Otsu method (8) or Kapur F kapur method (14) as the objective function then comparing the fitness value of X i and Q i and saving the global best solution with the highest fitness. We define the position of each agent based on the fitness value. In addition, we determine three probabilities of selection P 1 , P 2 , and P 3 using (21), (22), and (23) sequentially, and then, based on the probabilities, we calculate the position of each agent within the heap using (24). The agent's position (X) is updated using important D − ary heap-based operations, such as Heapify_Up(i), which is used to search for the superior node in the heap, and we insert the node at its correct position to preserve the heap characteristics, as demonstrated in Algorithm 1. Then, each agent upgrades its location frequently according to the best fitness value, and seeks the global optimum, as depicted in Algorithm 3. Optimization scenarios of implementing the proposed IHBO algorithm illustrated in Figure 1.

Segmentation phase.
In this phase, we generate the segmented image with the optimal threshold value in an image after setting x heap [1].value as the threshold value of the image. The pseudo-code of the proposed IHBO algorithm is illustrated in in Algorithm 4.

Performance evaluation of the proposed IHBO algorithm
Parameter settings. This section provides the estimation of the proposed IHBO algorithm. As we all know, adjusting parameters will certainly affect the performance of an algorithm. However, according to the suggestion of Arcuri et al. 69 , when comparing algorithm performance, all algorithm parameters should be kept at their default values, which come from their original papers, to ensure they are in a relatively optimal state. Moreover, we reduce the risk of better parametrization bias as each algorithm is set to its default values. Therefore, in this work, all algorithm parameters are kept at their default values. Thus, the performance of the proposed IHBO algorithm is evaluated over the IEEE Congress on Evolutionary Computation (CEC'2020) 70 as test problems. The CEC'2020 benchmark functions is utilized to test the performance of IHBO algorithm. Initially, this benchmark functions contained 10 test functions referred to as f 1 -f 10 . Consequently, function 1 is unimodal functions, functions 2-4 are multimodal functions, functions 5-7 are hybrid functions, and functions 8-10 are composition functions. Table 1 illustrates the parameters setting and mathematical formulation of the CEC'2020 benchmark functions; 'Fi*' refers to the best global value. Figure 2 illustrates a 2D visualization of the CEC'2020 benchmark functions to understand the differences and the nature of each problem.
Statistical results analysis of CEC'2020 benchmark. This  www.nature.com/scientificreports/ metrics. The standard deviation (STD) and mean of optimal solutions acquired by the proposed algorithm and all another algorithms utilized in the comparison is calculated. Furthermore, the qualitative metrics consists of average fitness history, convergence curve, and search history is used to evaluate the performance of the IHBO on the CEC'2020 test suite against seven well-known metaheuristic algorithms including the original HBO algorithm, SSA, MFO, GWO, SCA, HS, and EMO. Table 2 shows the STD and mean of the optimal value obtained from the proposed algorithm and the other competing algorithms for each CEC'2020 benchmark functions with 20 dimensional, and the optimal results of the STD and mean is minimum values in results. The results of the mean and STD of the proposed algorithm are proved superiority in solving seven CEC'2020 benchmark functions against to other competing algorithms.    Fig. 3. Boxplot is most important plots to describe data distributions into quartiles. This quartiles are the median of the first half of the data is first quartile, the second quartile is the median, the third quartile is median of the second half of the data, and the largest observation. The region among the first and third quartile is called the interquartile range and used to give an indication of spread in the data. The ends of the rectangles determine the lower and upper quartiles and a narrow boxplot refers to highest agreement among data. Figure 3 shows the boxplots of the proposed IHBO algorithm and illustrates the results of ten functions boxplot for Dim = 20. In reality, the results of proposed algorithm are proved superiority than all other competing algorithms on most of the test functions, but the performance of proposed algorithm is limited on F2, and F7.

Convergence curves analysis.
This subsection explains the convergence plots of the proposed algorithm with other competitor algorithms. Figure 4 illustrates the convergence plots of IHBO, HBO, SSA, GWO, MFO, HS, and SCA for the CEC'2020 benchmark functions. Furthermore, the proposed algorithm achieved optimal solutions and reached a stable point for most functions. Thus, IHBO can solve problems that require fast computation, such as online optimization problems. The proposed algorithm exhibited stable behavior, and its solutions converged easily in most of the problems it was tested on. Due to space limitations.  Figure 5 illustrates the qualitative analysis of the proposed IHBO algorithm. The first column illustrates a set of the CEC'2020 benchmark functions as plots in two-dimensional space. The second column illustrates the search history of search agents, from the first to the last iteration and display their exploitation behavior to realize the desired outcomes. The third column shows the average fitness history over 350 iterations, explaining the general behavior of the agents and the role that they play in the search of the best solution. According to average fitness history, all the history curves are decreasing, which means that the population improves at each iteration. The fourth column shows convergence curve and optimization history revealed the progress of fitness over a number of iterations. Optimization history is decreasing indicates that the solutions are optimized during iterations until reach the best solution.

Environmental and experimental requirements
This section presents the test images used for the experiments, then describes the empirical setup, and analyzes the results.   Table 3 displays the set of test images used.
Environmental setup. In this study, the proposed IHBO is compared with seven well-known metaheuristic algorithms including the original HBO, SSA, MFO, GWO, SCA, HS, and EMO. All competitor algorithms were applied and executed in Matlab 2015, and implemented on PC with 6G RAM running in a Windows 8.1 64-bit environment with an Intel Core I5 processor. The counterparts were executed 30 times per test image, number of iterations was set to 350, and population size is 50. The parameters settings of each algorithm were determined according to standard criteria and information found in previous literature (default values). The number of thresholds used was th 2 , th 3 , th 4 , and th 5 according to related literature 73 . The parameters settings of the IHBO and their values are presented in Table 4.
Evaluation metrics. Two metrics were utilized to estimate the performance of the IHBO algorithm. The first metric was used to evaluate the quality of the image, while the second metric was used to compare the edges of the segmented image. These metrics are important for evaluating the performance of the IHBO approach based on the Otsu and Kapur methods as objective functions. Statistical tests, such as the standard deviation (STD), Wilcoxon rank test, and average, were used to analyze the fitness of the proposed algorithm. We used the SSIM 48 , FSIM 74 , and PSNR 75 to measure the quality and stability of the image.
Structural similarity index (SSIM). The SSIM 48 index is a metric that is used to analyze the internal structures in a segmented image. A higher SSIM value indicates better segmentation of the original image due to the fact that structures in the two images match. The equation below describes the SSIM: The mean of the intensities of the original image I and segmented image Seg are µ I and µ Seg , respectively, and σ I and σ Seg are the standard deviations of the original image I and segmented image Seg, respectively. σ I,Seg is the covariance of the original image I and segmented image Seg, and c 1 and c 2 are two constants. The entire domain of the image is ω: Table 3. Set of test images. www.nature.com/scientificreports/ and G is the image's gradient magnitude and can be computed as follows: The vector's magnitude in v on n is E(v), and the local amplitude of scale n is A n (v) . The small positive number Peak signal-to-noise ratio (PSNR). The PSNR 75 is another metric used to evaluate the quality of segmentation by determining the difference between the quality of the original image and that of the segmented image. The PSNR is used to compare the original and segmented image using the root mean square error (RMSE) of each pixel, as expressed in (37). The PSNR can be defined as follows: where In (37), I and Seg are the segmented and original images of size M × N , respectively. A higher PSNR value indicates that there is higher similarity between the segmented and original images, which reflects a more effective segmentation process.

Experimental results and discussion
The experimental results are discussed in this section to evaluate the efficiency of the proposed algorithm.
Otsu results analysis. This subsection analyzes the outcomes of the IHBO based on the between-class variance as the fitness function, as proposed by Otsu. Table 7 illustrates the best threshold values obtained by applying the IHBO with the Otsu entropy as the objective function (8). Tables 5 and 6 present a graphical analysis of the thresholds, illustrating the resulting images of the IHBO with a different number of thresholds. Table 8 shows the computational time values of comparison algorithms obtained by Otsu's method. The IHBO proved its superiority in computational time compared to other competitive algorithms with 23 cases in 40 experiments and came in the first place. GWO came in second place with 10 experiments, while HBO come in third place with nine experiments, followed by EMO with two experiments. Finally, the MFO came in fifth place with only one experiment, and the remaining algorithms could not obtain the best computational time in any of the experiments. Table 9 illustrates the Otsu STD and average of the fitness results for the benchmark images. The IHBO demonstrated superiority in MTH by obtaining an optimal fitness values for 23 cases in 40 experiments. The HBO algorithm obtained the best fitness value in eight experiments, while the SCA obtained the optimal fitness value in five experiments and SSA come in fourth place with four experiments followed by MFO with three experiments. Finally, HS obtained the optimal fitness value in only one experiments and the remaining algorithms could not obtain the optimal fitness value in any of the experiments. Table 9 illustrates the STD values calculated for the 40 independent outcomes for each tested image with various thresholds. A lower STD value indicates that the algorithm is more stable. Table 10 presents the STD and mean PSNR for the benchmark images using the eight MAs. The IHBO was in first place in terms of the mean values of PSNR in 22 experiments. The SSA was in second place in seven experiments, while HBO was in third place, as it was superior in only six experiments. In fourth place was SCA with the best PSNR in only five experiments followed by MFO and HS with three experiments. Finally, the worst results were obtained by EMO which did not obtain the optimal values of the PSNR in any of the experiments. With respect to the STD, the IHBO was not the best alternative for lower dimensions (2 or 3 th). This is because the STD value was higher, which represents higher instability of the algorithm. However, MFO was a more unstable algorithm in terms of the PSNR. For the remaining approaches, the STD values followed the same tendency: lower for small dimensions and higher for four thresholds. However, the SSA was the least unstable algorithm, while the SCA was in second place. HBO was in the third place, HS was in fourth place, and the IHBO was in fifth place. Furthermore, GWO was in sixth place, and EMO was in seventh place. Table 11 illustrates the STD and mean of the FSIM obtained from 40 experiments. The results of the FSIM indicate that the IHBO obtained the highest FSIM in 22 experiments and was in first place, while the HBO was in second place in ten experiments. However, the SCA was in third place in eight experiments. SSA was in fourth place in two experiments, followed by HS and EMO, which appeared in fifth place in only one experiment. Finally, GWO and MFO came in last place in the experiments. The SCA was thus the best approach in terms of the STD because its values were lower in most experiments. The SSA came in second place, followed by EMO MxN . www.nature.com/scientificreports/ in third place. Then, the IHBO appeared in fourth place. GWO was in fifth place, followed by HBO. Finally, the least stable approaches were MFO and HS due to their high STD values in most cases. Table 12   www.nature.com/scientificreports/ SCA was the best method. In second place was the HBO, followed by EMO, which was in third place. The IHBO was in fourth place, while GWO was in fifth. Finally, MFO, HS, and SSA had no minimum STD values in the experiments. Table 7 illustrates the thresholds that were applied on the selected benchmark images. In Tables 5 and 6, the histograms are illustrated with the respective threshold values and the segmented images of the selected images using 2, 3, 4, and 5 thresholds. These results indicate that for some images, there was improvement in the quality of their contrast as the number of thresholds increased, particularly for the images Butterfly, Living Room, www.nature.com/scientificreports/ Jetplane, Lena, Pirate, Cameraman, Lake, and Bridge, presenting a higher amount of information in the image with the largest number of thresholds when compared with an image with only two thresholds. The most difficult histograms to segment were for Test 6, 9, and 10, relating to Bridge, Butterfly, and Barbara, respectively. The complexity was due to different numbers of pixels in the images, which could produce several classes or even make it impossible to select the optimal thresholds. Table 13 presents the p-values resulting from the Wilcoxon test for fitness using the Otsu fitness function. This table presents the difference between the proposed algorithm and the compared algorithms (HBO, SSA, MFO, GWO, SCA, HS, and EMO).
A difference between the SCA and MFO in comparison to the IHBO can be observed, which indicates that the proposed algorithm has a significant development. However, for the number of thresholds (nTh) = 5, the differences between the IHBO and most of the competing algorithms are clear by performing the comparison over 30 runs in each experiment. In the results, NaN indicates that the dataset to be compared is the same. This signifies that the algorithms obtained the same solution; thus, their results from the Wilcoxon test reveal that they are similar and that there are no differences between the methods.  Test 1   2  81 124  80 125  79 125  79 125  79 125  79 125  79 125  79 125   3  72 114 144  71 112 143  75 115 139  71 111 142  74 119 144  70 116 144  70 116 140  71 106 140   4  72 108 135 159  71 108 131 155  70 102 124 149  71 105 130 153  73 105 139 159  70 106 134 157  73 102 133 148  70 100 121 144   5  62 88 114 139 161 62 88 125 139 167 36 82 102 139 163 62 87 125 140 169 52 77 117 141 163 58 51 116 138 162 63 81 118 136 160 62 88 122 136 165   Test 2   2  65 122  65 122  65 122  63 122  65 122  65 121  65 122  65 122   3  50 88 128  65 122 129  55 89 126  58 88 124  59 92 129  53 85 124  59 91 127  50 86 123   4  48 85 118 150  48 85 118 150  50 88 118 151  48 82 115 149  46 81 117 144  42 87 110 141  45 88 115 150  39 80 112 140   5 48 81 107 133 160 48 85 102 145 159 48 85 108 140 159 48 88 101 143 160 48 83 106 138 157 44 83 101 138 155 48 85 103 145 150 48 79 102 141 Table 16, and were obtained by the IHBO using a fitness function such as Kapur entropy (14). Tables 14, and 15, present the histogram distribution of the benchmark images and segmented images with different numbers of thresholds produced by the IHBO. The results in Table 18 illustrate that the proposed algorithm with the Kapur entropy method proved outperform other algorithms in terms of SSIM (Table 18); in addition to, it outperformed other algorithms in terms of the mean FSIM (Table 20), PSNR (Table 19), and mean fitness. The values of the computational time of comparison algorithms obtained by Otsu's method are presented in Table 17. The IHBO came in first place with 24 cases in 40 experiments and proved its superiority in computational time compared to other competitive algorithms. HBO came in second place with 13 experiments, while GWO came in third place with ten experiments, followed by SSA with four experiments. Finally, the SCA came in fifth place with only one experiment, and the remaining algorithms could not obtain the best computational time in any of the experiments. Table 18 presents the STD and average fitness results of the Kapur method on the benchmark images. The IHBO was in first place by obtaining optimal fitness values with 24 cases in 40 experiments. The SCA was in second place in seven experiments, while the HBO was in third place in five experiments. SSA was in fourth place in three experiments, and HS was in fifth place in two experiments followed by GWO in sixth place in one experiments. Finally, EMO and MFO could not produce optimal fitness values. www.nature.com/scientificreports/ number of minimum STD cases, followed by IHBO in second place, HS in third place, MFO and the HBO in fourth place, SSA in fifth place, and SCA in sixth place. GWO had no optimal STD values. Table 20 provides the STD and mean FSIM. It can be observed that the IHBO was in first place in 17 experiments, while SCA was in second place in nine experiments. HS came in third place with six experiments, while MFO, HBO and GWO came in fourth place with three experiments, followed by SSA in fifth place with one experiment. Finally, EMO came in last place with no experiments. In terms of the STD, MFO came in first place with the maximum number of minimum STD cases, followed by SSA in second place, EMO and the IHBO in third place, HS in fourth place, HBO in fifth place, and SCA in sixth place. GWO had no optimal STD values.
The mean and STD of SSIM are presented in Table 21. The results indicate that IHBO was in first place in 19 experiments in terms of the SSIM, followed by SCA, which were in second place with seven experiments. The GWO and MFO were in third place in six experiments, while HBO and HS were in fourth place with three experiments, followed by SSA in fifth place with only one experiment. Lastly, EMO had no optimal experiments in terms of SSIM. According to the STD values, MFO came in first place with the maximum number of minimum STD cases. GWO, EMO, and SCA were in second place, while the IHBO was in third place. SSA and HBO were in fourth place, while HS was in fifth place.
Finally, Table 22 presents the p-values resulting from the Wilcoxon test for fitness using the Kapur fitness function. This table presents the difference between the proposed algorithm and the compared algorithms (HBO, SSA, MFO, GWO, SCA, HS, and EMO). The results in Table 22 indicate that the IHBO was different from the SCA and EMO but similar to the remaining algorithms. The exceptions occurred for nTH = 5, where in some cases the values exhibited differences as well as similarities (NaN values).
Human participants or animals. This article does not contain any studies with human participants or animals performed by any of the authors.

Conclusions and future works
Image segmentation is the most substantial pivotal phase that should be performed for image analysis and understanding. To handle this growing challenge, different methods using MTH, including feature-based, thresholdbased, and region-based segmentation, have been implemented. The most common technique used to perform and analyze image segmentation is threshold-based segmentation. This paper presented an improved variant of the Heap-based optimizer (HBO) called IHBO. The effectiveness of the proposed IHBO was estimated using the functions in the CEC'2020 benchmark functions, however, the proposed algorithm superiority on the competing algorithms regarding various statistical metrics. In addition, IHBO was applied to image segmentation using objective functions such as the Otsu and Kapur methods. The main target of IHBO is to determine the best thresholds that maximize the Otsu and Kapur methods. The IHBO was implemented on a set of test images with different characteristics, and the results were compared against seven well-known metaheuristic algorithms including the original HBO algorithm, SSA, MFO, GWO, SCA, HS, and EMO. The experimental results revealed that the IHBO algorithm outperformed all counterparts in terms of FSIM, SSIM, and PSNR. It should be noted that the IHBO results using the Otsu method provided better class variance in most metrics. However, when applying the Kapur method, the IHBO produced SSIM, FSIM, PSNR, and fitness values were better than those of all counterparts. The IHBO produced promising results because it preserved an effective balance between exploration and exploitation, and had the ability to avoid being trapped in local optima. For future work, there are many research directions in this field, such as studying the performance of the IHBO algorithm on different datasets, and other real-world complex problems. In addition, future work can study the hybridization of the original HBO with other metaheuristic or machine learning algorithms to automate the search process for the optimal number of thresholds in a specific image.       www.nature.com/scientificreports/  Test 1   2  95 164  96 163  96 163  96 163  96 164  96 163  96 163  97 162   3  24 96 164  23 96 163  23 96 163  23 96 163  23 109 161  23 96 163  23 96 163  23 95 164   4  23 80 125 173  23 80 125 173  23 80 125 173  23 80 125 173  23 74 130 184  23 80 125 173  23 80 125 173  47 48 139 140   5  23 62 94 135 177 23 77 119 159 190 23 62 94 135 177 23 62 94 135 177 21 76 113 135 166 23 71 109 144 180 23 71 109 144 180  23 111 115 119  138   Test 2   2  66 143  66 143  66 143  66 143  67 140  66 143  66 143  66 143   3  61 111 161  62 112 162  62 112 162  62 112 162  58 108 159  62 112 162  62 112